R Notebook

library(DT)

Attaching package: ‘DT’

The following object is masked from ‘package:Seurat’:

    JS

rr wolfrum <- readRDS(‘/projects/timshel/sc-scheele_lab_adipose_fluidigm_c1/data-wolfrum/wolfrum.compute.seurat_obj.rds’)

Which clusters are adipocytes/preadipocytes?

markers <- read.table('output/markergenes/wolfrum/markers_wolfrum.compute.seurat_obj.rds_seurat_clusters_negbinom', sep='\t', header=T)

Top 10 positive markers per cluster

Which clusters are the preadipocytes? Check by plotting some of the DE genes between T1T2T3 and T4T5 in the 180831 dataset.

top10_T1T2T3 <- T1T2T3[order(-T1T2T3$avg_logFC),][1:10]

Error in `[.data.frame`(T1T2T3[order(-T1T2T3$avg_logFC), ], 1:10) : 
  undefined columns selected

How do these genes look in the 180831 data?

How are they expressed in the Wolfrum data?

Not all genes are expressed clearly in specific clusters. It looks like cluster 22, 21, 5, 14, 23, 11 and 10 are the preadipocytes.

Also plot the top 10 ECM and Metabolic markers.

markers_u_l <- markers_u_l[order(-markers_u_l$avg_logFC)]

Error in `[.data.frame`(markers_u_l, order(-markers_u_l$avg_logFC)) : 
  undefined columns selected

How do these genes look in the 180831 data?

And how are they expressed in the Wolfrum data?

plots <- FeaturePlot(wolfrum, features=c(as.vector(markers_u$gene)[1:10], as.vector(markers_l$gene)[1:10]), pt.size=1, combine=F)

The following requested variables were not found: RP11-572C15.6

plot_grid(plotlist=plots, ncol=2)

The U and L branch markers are clearly expressed in the two clusters in the top right of the UMAP plot.

Check how many of the marker genes are found in each cluster.

kable(df[order(-df$total_overlap),]) %>%
  kable_styling(bootstrap_options = "striped", full_width = F)

<table class="table table-striped" style="width: auto !important; margin-left: auto; margin-right: auto;">
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> total_overlap </th>
   <th style="text-align:right;"> overlap_P </th>
   <th style="text-align:right;"> overlap_L </th>
   <th style="text-align:right;"> overlap_U </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> 11 </td>
   <td style="text-align:right;"> 27 </td>
   <td style="text-align:right;"> 13 </td>
   <td style="text-align:right;"> 18 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 14 </td>
   <td style="text-align:right;"> 23 </td>
   <td style="text-align:right;"> 10 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 14 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 22 </td>
   <td style="text-align:right;"> 22 </td>
   <td style="text-align:right;"> 15 </td>
   <td style="text-align:right;"> 7 </td>
   <td style="text-align:right;"> 5 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 23 </td>
   <td style="text-align:right;"> 20 </td>
   <td style="text-align:right;"> 9 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 11 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 10 </td>
   <td style="text-align:right;"> 18 </td>
   <td style="text-align:right;"> 12 </td>
   <td style="text-align:right;"> 9 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 5 </td>
   <td style="text-align:right;"> 10 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 8 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 24 </td>
   <td style="text-align:right;"> 10 </td>
   <td style="text-align:right;"> 8 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 2 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 7 </td>
   <td style="text-align:right;"> 6 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 13 </td>
   <td style="text-align:right;"> 6 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 2 </td>
   <td style="text-align:right;"> 5 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 3 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 15 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 19 </td>
   <td style="text-align:right;"> 4 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 18 </td>
   <td style="text-align:right;"> 3 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 4 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 6 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 8 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 12 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 1 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 20 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 21 </td>
   <td style="text-align:right;"> 2 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 1 </td>
   <td style="text-align:right;"> 1 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 1 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 9 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 16 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 17 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 25 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> 26 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
   <td style="text-align:right;"> 0 </td>
  </tr>
</tbody>
</table>

Table shows the percentage of genes with avgLogFC > 0.7 that was found in the cluster. Sort on columns to get the top clusters.

logfc_0.5_genes_100 <- get_gene_overlap_per_cluster(0.5, 100)
datatable(logfc_0.5_genes_100)

Check the top 5 clusters per branch for different n_genes and min logFC

logfc <- c(0.25, 0.5, 0.7)
n_genes <- c(100, 100, 100)
for (i in 1:length(logfc)){
  print(paste('min_logFC: ', logfc[i], ' | n_genes: ', n_genes[i], sep=''))
  df <- get_gene_overlap_per_cluster(logfc[i], n_genes[i])
  print_top_clusters(df)
}

[1] "min_logFC: 0.25 | n_genes: 100"
[1] "P: 24, 11, 22"
[1] "L: 11, 10, 23"
[1] "U: 14, 24, 23"
[1] "min_logFC: 0.5 | n_genes: 100"
[1] "P: 24, 11, 22"
[1] "L: 11, 10, 23"
[1] "U: 14, 23, 24"
[1] "min_logFC: 0.7 | n_genes: 100"
[1] "P: 24, 11, 19"
[1] "L: 11, 10, 23"
[1] "U: 14, 23, 22"

There is some overlap between branches which is good. Cluster 11 is shared between P and L, cluster 14 is shared between P and U and cluster 23 is shared between L and U. This would indicate that cluster 24 contains most immature preadipocytes, cluster 10 contains most mature L branch cells and cluster 22 contains most mature U branch cells.\

Based on the results:\ P = 24\ L = 11\ U = 14\

Hypothesis: cluster 23 represents preadipocytes at the start of differentation (the cell states between T3 and T4 in 180831 data that we missed). Cluster 5 represents even more mature metabolic cells and cluster 11 represents more mature ECM cells.\

Clusters 22 shares most genes with the P branch. Cluster 23 most with the U branch. (see datatable above). These could also represent the preadipocytes at start of differentiation or the cells that transfer back to progenitor cells. \

Seurat integration of Wolfrum and 10x-180831 data

These results also confirm that the L branch is closest to cluster 11 and U is closest to the U branch.\

Seurat integration of Wolfrum preadiopcyte subset and 10x-180831 data

integrated.subset@meta.data[which(is.na(integrated.subset@meta.data$branch)), 'dataset'] <- 'Wolfrum'

Error in integrated.subset@meta.data[which(is.na(integrated.subset@meta.data$branch)),  : 
  object 'integrated.subset' not found

Predict cell types with Seurat’s TransferData

wolfrum.predicted_labels <- readRDS('output/seurat_objects/wolfrum/wolfrum.predicted_labels_180831.rds')

Used pca, pca.project and cca as dimred for FindTransferAnchors and IntegrateData. \

Scores for P cells

Scores for U cells

Scores for L cells

For all predictions, change predicted id to NA if max score is below a certain threshold.

assign_labels <- function(colname, threshold=0.5){
  pred_ids <- unlist(as.vector(apply(wolfrum.predicted_labels@meta.data[,c(paste(colname,'.prediction.score.max', sep=''), paste(colname, '.predicted.id', sep=''))], 1, function(x){
    if (x[[1]] < threshold){
      return(NA)
    } else{
      return(x[[2]])
    }
  })))
  return(pred_ids)
}

for (col in c('predictions_pca_project', 'predictions_pca', 'predictions_cca.')){
  for (t in c(0.5, 0.7, 0.9, 0.95, 0.99)){
    preds <- assign_labels(col, t)
    wolfrum.predicted_labels <- AddMetaData(wolfrum.predicted_labels, preds, col.name=paste(col, 'predicted_label', t, sep='.'))
  }
}

Threshold for prediction = 0.5

plot_grid(
  UMAPPlot(wolfrum.predicted_labels, group.by='predictions_pca.predicted_label.0.5'),
  UMAPPlot(wolfrum.predicted_labels, group.by='predictions_pca_project.predicted_label.0.5'),
  UMAPPlot(wolfrum.predicted_labels, group.by='predictions_cca..predicted_label.0.5'), ncol=2
)

Threshold for prediction = 0.7

Threshold for prediction = 0.9

Only show preadipocytes

Figures

labels <- unlist(lapply(wolfrum.predicted_labels@meta.data$predictions_pca.predicted_label.0.7, function(x){
  if (is.na(x)){
    return('Non-maching cells')
  } if(x == 'ECM'){

Error: unexpected 'if' in:
"    return('Non-maching cells')
  } if"

save_plot("figures/figures_paper/main_figures/Figure_wolfrum/UMAP_wolfrum_predicted-labels_180831.pdf", p_predictions, base_width=8, base_height=5)

Error in grDevices::pdf(file = filename, ..., version = version) : 
  cannot open file 'figures/figures_paper/main_figures/Figure_wolfrum/UMAP_wolfrum_predicted-labels_180831.pdf'

p <- UMAPPlot(wolfrum.predicted_labels, group.by='RNA_snn_res.0.8', label=T, no.axes=T, no.legend=T) + theme(legend.position = "none", axis.text = element_blank(), axis.ticks = element_blank(), axis.title = element_blank(), plot.margin=grid::unit(c(0,0,0,0), "mm"))

The following functions and any applicable methods accept the dots: CombinePlots

#save_plot("figures/figures_paper/main_figures/Figure_wolfrum/UMAP_wolfrum_clusters.pdf", p_clusters, base_width=6, base_height=5)

#featureplots
#widht=13
#ucp2_fp <- featureplots_leg$UCP2 + scale_color_gradient(name='Expression', low='gray', #high='blue', guide='colorbar', limits=c(0,5)) + theme(plot.title=element_blank(), #legend.title=element_text(size=20), legend.text=element_text(size=20), legend.key.height = #unit(1.3, 'cm'))


#dcn_fp <- featureplots_noleg$DCN + theme(plot.title=element_blank())

adipoq <- FeaturePlot(wolfrum.predicted_labels, features='ADIPOQ') + NoLegend() + NoAxes() + theme(plot.title = element_text(size=20)) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5))

lipe <- FeaturePlot(wolfrum.predicted_labels, features='LIPE') + NoLegend() + NoAxes() + theme(plot.title = element_text(size=20)) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5))

apod <- FeaturePlot(wolfrum.predicted_labels, features='APOD') + NoLegend() + NoAxes() + theme(plot.title = element_text(size=20)) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5))

dcn <- FeaturePlot(wolfrum.predicted_labels, features='DCN') + NoAxes() + theme(plot.title = element_text(size=20), legend.text=element_text(size=20), legend.key.height = unit(1.3, 'cm')) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5)) 

#save_plot("../figures/figures_paper/main_figures/Figure_wolfrum/featureplots.pdf", g, base_width=16, base_height=4)

EBF2 and LEP

pparg <- FeaturePlot(wolfrum.predicted_labels, features='PPARG') + NoAxes() + theme(plot.title = element_text(size=20), legend.text=element_text(size=20), legend.key.height = unit(1.3, 'cm')) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar') 

Scale for 'colour' is already present. Adding another scale for
'colour', which will replace the existing scale.

save_plot("figures/figures_paper/main_figures/Figure_wolfrum/UMAP_adipocyte_clusters.pdf", p, base_width=6, base_height=5)

Error in grDevices::pdf(file = filename, ..., version = version) : 
  cannot open file 'figures/figures_paper/main_figures/Figure_wolfrum/UMAP_adipocyte_clusters.pdf'

ebf2 <- FeaturePlot(wolfrum.predicted_labels, features='EBF2', pt.size=0.5) + NoLegend() + NoAxes() + theme(plot.title = element_text(size=20)) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5))

Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the
existing scale.

pparg <- FeaturePlot(wolfrum.predicted_labels, features='PPARG', pt.size=0.5) + NoAxes() + theme(plot.title = element_text(size=20), legend.text=element_text(size=15), legend.key.height = unit(1.3, 'cm')) + scale_color_gradient(name='Expression', low='gray', high='blue', guide='colorbar', limits=c(0,5)) 

Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the
existing scale.

g <- plot_grid(
  ebf2, pparg, ncol=2, rel_widths = c(1, 1.3)
)
g

plots <- VlnPlot(wolfrum.predicted_labels, features=c('ADIPOQ', 'LIPE', 'PLIN4', 'FABP4', 'ADIRF', 'APOD', 'MGP', 'DCN', 'CCDC80', 'PLAC9'), group.by='RNA_snn_res.0.8', pt.size=-1, combine=F)

for (i in 1:length(plots)){
  if (i == length(plots)){
    plots[[i]] <- plots[[i]] + NoLegend() + 
      theme(plot.title=element_blank(), 
            axis.title.y=element_blank(), 
            axis.line.x=element_blank(),
            axis.text.x=element_text(angle=0, size=12),
            plot.margin = unit(c(0, 0, 0, 0), "cm")) + 
      labs(x='Cluster')
  } else {
    plots[[i]] <- plots[[i]] + NoLegend() + 
      theme(plot.title=element_blank(), 
            axis.title.y=element_blank(), 
            axis.line.x=element_blank(),
            axis.ticks.x=element_blank(),
            axis.text.x=element_blank(),
            axis.title.x=element_blank(),
            plot.margin = unit(c(0, 0, 0, 0), "cm"))
  }
}

vlnplts <- plot_grid(plotlist=plots, ncol=1, rel_heights=c(1,1,1,1,1,1,1,1,1,1.6))
vlnplts

#save_plot("figures/figures_paper/main_figures/Figure_wolfrum/violinplots.pdf", p, base_width=6, base_height=8)

Supplementary figure

integrated <- readRDS('output/seurat_objects/wolfrum/wolfrum.180831.integrated.rds')
sfig <- plot_grid(
  UMAPPlot(integrated, group.by='dataset'),
  UMAPPlot(integrated, group.by='State.labels', cols=colormap.branches),
  UMAPPlot(integrated, group.by='RNA_snn_res.0.8', label=T) + NoLegend(), ncol=2
)
sfig

#save_plot("figures/figures_paper/supplementary_figures/wolfrum/integration.wolfrum.180831.pdf", sfig, base_width=12, base_height=8)